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We compare simulations performed using the Nose-Hoover and the Langevin thermostats 
with the Hamiltonian dynamics of a long-range interacting system in contact with a reser- 
voir. We find that while the statistical mechanics equilibrium properties of the system are 
recovered by all the different methods, the Nose-Hoover and the Langevin thermostats 
fail in reproducing the nonequilibrium behavior of such Hamiltonian. 
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1. Introduction 

Long-range characterizes the interactions of a number of different physical sys- 
tems like, e.g;, plasmas, wave-matter systems, gravitational systems, Bose-Einstein 
condensates In these cases, the customary assumptions of statistical mechanics 
are put into question because of the inapplicability of the Boltzmann transport 
equation.^! In fact, nonequivalence between the microcanonical and the canonical 
ensemble approaches,^ nonergodicity and topological nonconnetctivitjB has been 
detected for long-range system. Classical long-range Hamiltonian models assume 
then a central role in order to compare the dynamical behavior of macroscopic 
phase functions like the system's energy, its temperature, or its magnetization, with 
the correspondent predictions of statistical mechanics. It is known that the dynam- 
ics of long-range Hamiltonians displays long-living quasi-stationary states (QSS) in 
microcanonical (fxC) simulations, i.e., when the system is isolated. This aspect has 
been studied in details by several groups in the last decadeP On the other hand, at 
least in terrestrial-scale experiments, the system cannot be considered isolated. It 
is then interesting to see if QSSs are reproduced in more "experimental" settings, 
especially in view of some theoretical results, based on the Langevin equation, that 
seem to rule out such a possibility.^ 

In Ref.[?l[H]we addressed this issue by introducing a Hamiltonian setup in which 
the long-range system is coupled with a thermal reservoir through microscopic in- 
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teractions. We discovered the persistence of long-lasting QSSs whose life-time de- 
pend on the system size and on the coupling strength between the system and the 
reservoir In this Paper we further investigate this point by comparing standard 
methods for the simulation of a thermal bath interacting with the system, namely 
the Nose-Hoover (NH) and the Langevin (LA) thermostats, with the above Hamil- 
tonian. 

An important observation with respect to the NH and the LA schemes is that 
both algorithms implicitly assume some equilibrium features. In the first case, a 
single degree of freedom with "effective mass" Q is added to the system in order 
to simulate a thermal bath. This additional degree of freedom has the capability 
of adsorbing and releasing an arbitrary amount of energy, and its temporal scale 
is appropriately redefined in such a way to generate the Boltzmann-Gibbs equi- 
librium canonical distribution for the system!^ The parameter Q can be used for 
tuning a better convergence of the algorithm. At difference, the Langevin approach 
is based on the assumption of a well defined separation between the time-scales of 
the dynamics of the diffusive particle (slow dynamics) and that of the underlying 
thermal bath (fast dynamics). Because of this, the bath is assumed to be in thermal 
equilibrium at all the integration steps and via the equipartition theorem the damp- 
ing and the stochastic coupling constant characterizing the diffusive behavior are 
related by a specific, temperature-dependent, fluctuation-dissipation relation.^ As 
a consequence of these implicit assumptions, it is not guaranteed that the NH and 
the LA integration strategies can be safely applied out of equilibrium, especially for 
a system which is known to produce unconventional effects in fiC simulations. 

In the following, we take an "empirical" attitude by implementing the NH and 
the LA integration schemes for the simulation of a long-range Hamiltonian system 
in contact with a thermal bath. We find that while the equilibrium behavior of 
the system is equivalently recovered by the different approaches, the NH and the 
LA thermostats do not properly account for the noncquilibrium features of the 
Hamiltonian dynamics. 

2. Equilibrium simulations 

The long-range interacting Hamiltonian considered in Refs.m[6j|7l|8]is called Hamil- 
tonian Mean Field (HMF) model and it can be thought as a set of M globally 
coupled Ay-spins with Hamiltonian 

M p l M 

Hhmf = Khmf + Vhmf = X^ + ^7^[ 1 ~ cos ( 01 ~ 9 ^ ' 0-) 

where 6i £ [0, 2tt) are the spin angles assumed with unit momentum of inertia 
and U £ M their angular momenta (velocities). Its equilibrium statistical me- 
chanics solution predicts an high-energy disordered phase separated from a low- 
energy ordered one by a second order transition occurring at the specific energy 
&hmf = Ehmf/M — 0.69 (we use dimensionless units). The order parameter is the 
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Fig. 1. Caloric curve. The solid line is the Boltzmann-Gibbs equilibrium solution and the dashed 
line is the prolongation of the ordered phase to subcritical energies. Empty symbols are the average 
value of enMF{t) a t equilibrium. Full circles correspond to the QSS studied in the paper and to 
the subsequent microcanonical and canonical equilibrium obtained as t — * oo. 

magnetization of the system nriHMF = | X)f=i( cos ^' sin0j)|/M, and the presence 
of the kinetic term endows the spin system with a proper Hamiltonian dynamics 
in which one can define the temperature T as twice the specific kinetic energy, 
Thmf = 2K H mf/M. Notice the relation chmf = (T H mf + 1 - m 2 HMF )/2. 

The Hamiltonian thermal bath (HTB) considered in Refs. [71 [5] is given by the 
full Hamiltonian system H = Hrmf + Htb + Hi, 

N j 2 N 

H tb= £ |+ £ [l-cosffc+i-ft)], (2) 

i=M+l i=M+l 
M S 

ffi = e££[l-cos(0i-0 r , (<) )], (3) 

i=l s=l 

where Htb and Hi are respectively the Hamiltonian of the thermal bath and that of 
the interaction between HMF model and thermal bath (JV 3> M). As the coupling 
constant e vanishes, the fiC dynamics of the HMF is recovered (see Ref.UJfor details). 

A similar circumstance is valid for the LA thermostat where the equation of 
motion for the HMF model are 

M 

9 i = -7e i -J2M6i-0 j ) + V2lf^(t), i 1.2 \l. (4) 

j=i 

In Eq. ([J} £i(t) is a Gaussian white noise characterized by a zero average = 

and correlation (£,i(t)£,j(t')) — 5ijS(t — t'). Indeed, in the limit 7 — > 0, Eqs. ^ reduce 
to the /xC Hamiltonian equations of the HMF model. 

In Fig. [1] we present the results of the simulations in the different setups at 
equilibrium. These simulations are obtained by setting initial conditions close to 
equilibrium for the HMF model. Specifically, we used a Maxwellian distribution of 
velocities and an initial value of tuhmf and of Thmf close to those of equilibrium 
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Fig. 2. QSSs in fiC simulations (a) and in HTB simulations with e = 0.01 (b). Dashed lines are 
the equilibrium values for m 2 HMF - See Refs. 7,8 for further details. 

for the given fixed energy (/LtC) or thermal bath temperature (HTB, NH and LA). 
After the relaxation to equilibrium, we verified that all the different dynamics give 
the same results for the values of the phase functions chm f > Thmf and tuhmf- 



3. Nonequilibrium simulations 

The nonequilibrium analysis is performed by changing the initial conditions for 
the dynamical variables of the HMF model, while keeping all the other parameters 
unchanged with respect to the equilibrium results. Consistently with other studies 
reporting the existence of QSSsp^ we have chosen at time t = 0a delta distribution 
for the angles {phmf(6, 0) = 5(0) so that m 2 HMF (G) = 1] and a uniform distribution 
for the velocities \phmf(1,0) = 1/2^, I € [— I, I], with I ~ 2.03]. In this way, the 
initial energy of the HMF model is set to the subcritical value chmf^) — 0.69. 

The /LtC simulations reveal in this case the existence of a violent relaxation 
process^ (for a time of order 1) followed by a QSS which can be displayed e.g. 
by plotting the time dependence of m 2 HMF (Fig. [2^). The QSS life-time diverges 
in the thermodynamic limit M — > oo and in this limit rn 2 HMF vanishes. The same 
kind of results are obtained using the HTB, although now the QSS life-time dimin- 
ishes as e increases (Fig. [2]o) . Notice that in the /j,C simulations the system relaxes 
to equilibrium at fixed energy, whereas in the HTB ones the relaxation is at fixed 
thermal bath temperature. This produces a consistent difference in the equilibrium 
values of m 2 HMF (see also Fig. [T]). 

Unlike the results in the previous section, a NH integration scheme implemented 
with nonequilibrium initial conditions for the HMF model does not always guarantee 
the final convergence to equilibrium (Fig. [3]) . We found that only when Q is larger 
then the value of the system size M the convergence to equilibrium is realized. Still, 
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Fig. 3. NH simulations with far-from-equilibrium initial conditions, (a) With Q = 1 the magne- 
tization does not converge to equilibrium, (b) With Q = M = 10 3 the magnetization converges to 
equilibrium without displaying QSSs. 



the nonequilibrium dynamics is characterized by fluctuations of the phase functions 
which display no relation with the Hamiltonian simulations. Another drawback of 
the NH method is that if initially the HMF model has vanishing total momentum, 
this quantity remains zero during all the integrations steps. 

The analysis of the LA simulations reveals some interesting new results. In this 
case, as for the HTB, the convergence to equilibrium is observed for any value of 
7 > and the total momentum of the HMF model fluctuates during the simulation, 
as it is expected. Also, the violent relaxation process is coherently reproduced by the 
LA simulations and a QSS follows for which m 2 HMF — > as M grows. Nonetheless, 
the QSS life-time appears to be independent from the system size M for any value 
of 7 > (Fig. [4^). This life-time also shows an interesting dependence on 7. While 
the violent-relaxation time is of order 1 independently on 7, the crossover time from 
the QSS to the equilibrium scales as jt (Fig.|3J}). This scaling law implies an infinite 
life-time of the QSS in the fiC limit 7 — > 0, independently on M. Since such a result 
is in contrast with purely Hamiltonian [iC simulations (Fig. [2^), it suggests the 
presence of a discontinuity in 7 = 0. 

It is interesting to recall that a stability analysis of the Fokker-Planck equation 
derived from Eq. ((4|) shows that anomalous, non-Maxwellian, velocity probability 
density functions are (neutrally) stable only in the fiC limit 7 = oPThe somehow 
unexpected^ presence of QSSs in LA simulations may be related to the fact that 
during the QSS the HMF model does not thermalize with the thermal bath (see 
Ref.HQlfor details). 

In conclusion, by showing a specific example in which the NH and the LA ther- 
mostats simulations do not agree with the correspondent fully Hamiltonian ones, 
our findings constitute a general warning against the straightforward application of 
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Fig. 4. LA simulations with far-from-equilibrium initial conditions, (a) As M increases, the 
QSS's life-time remains constant, (b) The crossover time between QSS and equilibrium scales as 
It. 

equilibrium-based algorithms for the description of the statistical nonequilibrium 
behavior. 
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